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ABSTRACT 
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This  report  describes  automatic  differentiation,  which  is  neither 
symbolic  nor  approximate,  for  single  functions  of  one  real  variable.  The 
rules  of  evaluation  and  differentiation  are  combined  into  an  ordered-pair 
arithmetic  similar  to  complex  arithmetic,  but  slightly  simpler.  Evaluation  of 
the  formula  for  a  function  in  this  arithmetic  yields  both  the  values  of  the 
function  and  its  derivative,  without  a  formula  for  the  derivative  of  the 
function,  and  without  numerical  approximation,  since  this  arithmetic  is  based 
on  the  well-known  rules  for  differentiation.  The  properties  of  this 
arithmetic  are  examined,  and  illustrated  by  simple  examples.  Subroutines  are 
given  for  differentiation  arithmetic  both  on  a  hand-held  programmable 
calculator,  and  in  the  microcomputer  language  Pascal-SC.  An  application  of 
this  arithmetic  to  the  solution  of  equations  by  Newton's  method  is  given, 
using  a  Pascal-SC  program. 
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SIGNIFICANCE  AND  EXPLANATION 


Most  people  think  that  evaluation  of  derivatives  on  a  computer  has  to  be 
either  symbolic  or  numerical,  that  is,  approximate.  By  combining  the  rules  of 
differentiation  with  the  rules  for  evaluation,  however,  it  is  possible  to 
perform  differentiation  using  only  the  formula  (or  subroutine)  for  the 
function  involved.  It  is  not  necessary  to  produce  a  formula  for  the 
derivative  of  the  function  to  obtain  the  exact  value  of  its  derivative  for  any 
given  values  of  its  arguments.  The  way  this  method  of  automatic 
differentiation  works  is  not  immediately  obvious,  because  everyone  is  taught 
in  school  that  the  first  step  in  evaluating  the  derivative  of  a  function  is  to 
derive  the  formula  for  the  derivative,  and  the  only  alternative  to  this  is  to 
approximate  the  derivative  by  a  difference  quotient.  However,  the  definition 
of  the  dervative  and  the  rules  derived  from  it  actually  define  values  of  the 
derivative,  not  formulas,  so  that  these  rules  can  be  used  to  define  what  is 
called  differentiation  arithmetic  here  for  single  real-valued  functions  of  one 
real  variable.  This  is  an  ordered-pair  arithmetic,  similar  to  complex 
arithmetic  but  actually  a  little  simpler.  When  the  formula  for  the  function 
is  evaluated  in  this  arithmetic,  the  value  of  the  function  and  its  derivative 
are  obtained  automatically,  without  symbolics  and  without  approximations. 
Differentiation  arithmetic  is  easy  to  program  for  hand-held  calculators  as 
well  as  computers,  and  some  examples  are  given. 

Automatic  differentiation  was  developed  first  in  the  general  setting  of 
functions  of  several  variables,  and  for  Taylor  series  expansions  of  functions 
of  single  variables.  The  special  case  considered  here  was  presented  to  the 
Wisconsin  Section  of  the  Mathematical  Association  of  America  in  the  hope  that 
this  simple  setting  would  help  to  explain  away  some  of  the  current 
misconceptions  about  differentiation  on  the  computer.  Automatic 
differentiation  is  as  accurate  as  symbolic  differentiation,  which  requires  a 
lot  of  software,  but  much  faster.  It  is  also  suitable  for  parallel 
computation.  Numerical  differentiation  gives  only  approximate  answers,  while 
automatic  differentiation  gives  correct  results  at  about  the  same  speed  or 
even  faster. 


The  responsibility  for  the  wording  and  views  expressed  in  this  descriptive 
summary  lies  with  MRC,  and  not  with  the  author  of  this  report. 


THE  ARITHMETIC  OP  DIFFERENTIATION 


L.  B.  Rail 


1.  A  simple  calculus  problem.  Traditional  ways  of  thinking  and  doing  things  exert  a 
strong  hold  in  mathematics,  as  well  as  in  other  human  endeavors.  For  example,  generations 
of  students  were  taught  to  solve  linear  systems  of  equations  by  using  determinants  and 
Cramer's  rule  [2].  It  was  not  until  people  actually  tried  solving  linear  systems  using 
desk  calculators  and  later  computers  that  more  efficient  methods  based  on  Gaussian 
elimination  came  into  general  use  [7],  (3],  Even  today,  it  is  disconcerting  to  encounter 
people  who  think  that  the  "only  way"  to  solve  linear  systems  is  by  the  Cramer's  rule  they 
learned  in  school. 

Calculus  is  another  subject  with  a  long,  rich  history  and  strong  traditions.  As  an 
example  of  a  traditional  approach  to  a  simple  problem  in  calculus,  suppose  you  are  given  a 
function  defined  by  a  formula,  say 


(1.1) 


F(X)  - 


X2  +  2X  -  3 
X  *  2 


(X  -  1).(X  +  3) 
X  *  2 


and  you  want  to  find  the  value  of  the  derivative  F'(X)  for  a  given  value  of  X.  Most  people 
have  been  taught  to  solve  this  problem  by  first  deriving  a  formula  for  the  derivative. 
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(1.2) 


F*  (X) 


X2  +  4X  +  7 
(X  +  2)^ 


+  (X  +  2)'1  ' 


and  then  using  the  usual  rules  of  arithmetic  to  evaluate  F'(X).  This  method,  like  Cramer's 
rule,  is  theoretically  correct,  and  not  complicated  for  this  example.  However,  this  is  not 
the  only  way  to  calculate  the  values  of  derivatives.  It  turns  out  that  it  is  possible  to 
obtain  both  the  values  F(X)  and  F'(X)  simply  by  evaluating  the  formula  for  F(X)  using  a 
different  kind  of  arithmetic,  called  differentiation  arithmetic.  A  formula  for  F'(X)  is 
not  required,  contrary  to  common  belief.  Differentiation  arithmetic  is  an  ardered-pair 
arithmetic,  similar  to  complex  arithmetic.  The  rules  for  this  arithmetic  will  be  defined 
lateri  first,  it  is  important  to  analyze  the  key  idea  of  function  evaluation  in  ordinary 
arithmetics,  and  see  what  this  process  has  in  common  with  differentiation. 

2.  Evaluation  of  functions.  The  usual  way  to  evaluate  a  rational  function,  (1.1)  for 
example,  is  to  break  it  down  into  a  sequence  of  additions,  subtractions,  multiplications, 
and  divisions,  and  then  apply  the  rules  for  arithmetic  to  each.  Given  a  simple  formula, 
this  is  often  done  instinctively,  without  ouch  thought.  If  a  hand-held  calculator  is  being 
used,  however,  it  is  useful  to  think  out  the  sequence  of  arithmetic  operations  being  used, 
and  perhaps  write  them  down,  especially  if  a  complicated  formula  is  being  evaluated.  Such 
a  sequence  is  called  a  code  list  for  F(X)  (10),  [12].  For  example,  evaluation  of  (1.1)  can 
be  done  in  the  following  steps: 


T1  X  -  1» 

T2  X  +  3; 

(2.1)  T3  T1 *T2; 

T4  X  +  2; 

F(X)  T3/T4| 


or  sosw  other  equivalent  sequence.  The  intermediate  quantities  T1,T2,T3,T4  are  each  the 
result  of  a  single  arithmetic  operation. 

Another  way  to  visualize  the  evaluation  of  FIX)  is  by  means  of  a  Kantorovich  graph 
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(5) ,  fll),  112),  shown  in  Figure  2.1. 

The  Kantorovich  graph  is  a  directed  graph  in  which  information  (numerical  values  in 
this  case)  flows  along  edges  from  top  to  bottom,  and  the  indicated  arithmetic  operations 
are  performed  at  the  nodes.  This  graph  shows  which  operations  somt  be  performed  before 
others  can  be  done,  or,  from  another  point  of  view,  the  operations  which  can  be  performed 
independently  or  in  parallel.  Various  code  lists  equivalent  to  (2.1)  can  be  read  off  the 
graph.  In  the  early  days  of  computing,  the  tedious  process  of  preparation  of  code  lists  to 
evaluate  functions  was  done  by  the  programmer.  Modern  computer  languages  now  contain 
formula  translators,  which  automatically  analyze  formulas  such  as  (1.1)  by  a  fixed 
algorithm  and  produce  sequences  of  instructions  similar  to  (2.1)  for  the  machine  to 
execute. 


Figure  2.1.  A  Kantorovich  Graph  for  the  Function  (1.1). 


Once  the  code  list  (2.1)  has  been  formed,  evaluation  of  P(X)  for  some  value  of  X,  say 
X  •  3,  is  a  trivial  matter;  3  is  substituted  for  X  and  the  indicated  sequence  of  operations 
are  performed: 


(2.2) 


T1  i-  3  -  1  -  2; 

T2  3  +  3  -  6; 

T3  2*6  -  12; 

T4  3  +  2  -  5; 

F(3)  12/5  -  2.4  . 


Of  course,  code  lists  are  not  the  only  vay  to  express  evaluation  of  a  function.  Another 
popular  method  is  to  express  F(X)  in  "reverse  Polish"  notation,  for  example,  by 

(2.3)  F(X)  Xt1-Xt3+*Xt2+/; 

where  t  represents  "entering"  the  preceding  value,  or  "pushing*  it  onto  a  stack. 

3.  Differentiation  of  functions.  Elementary  differential  calculus  is  based  on  the 
definition 

<3  1\  .  Um  F(X  +  AX)  -  F(X) 

'  1  AX  ♦  0  AX 

However,  this  definition  is  seldom  used  for  actual  differentiation  of  a  function  such  as 
(1.1).  Instead,  it  is  used  to  derive  roles  for  differentiation  of  sums,  differences, 
products,  quotients,  and  other  functions  encountered  in  mathematical  analysis.  Once  these 
rules  have  been  established,  differentiation  proceeds  in  much  the  same  way  as  function 
evaluation:  The  function  is  expressed  as  a  sequence  of  additions,  subtractions, 
multiplications,  divisions,  etc.,  and  the  appropriate  rules  are  applied.  This  technique  is 
valid  as  a  consequence  of  the  chain  role,  which  is  derived  from  the  definition  (3.1).  To 
see  how  this  works  for  our  example,  apply  the  rules  for  differentiation  to  the  code  list 
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<2-1 )  s 


Tl* 

!-  1> 

T2* 

1l 

(3.2) 

T3* 

s-  Tl  *T2 1  +  T2T1'» 

T4 1 

1j 

F 1  (X) 

(T4 *T3 '  -  T3T4M/T42. 

This, 

along  with  (2.2),  gives  all 

the  information  needed  to  evaluate  F'(3).  The  steps  are: 

Tl  • 

i-  1> 

T2  • 

1> 

(3.3) 

T3* 

2*1  +  6*1  -  8; 

T4 1 

1» 

FM3) 

(5*8  -  12*1 )/52  -  28/25  -  1.12  . 

• 

•H 

£ 

gives  the  correct  value  of  F 

"(3)  without  the  use  of  a  formula  for  F'(X).  This  process 

ia  called  antontlc  differentiation  [11],  [12]  in  order  to  distinguish  it  from  symbolic 
differentiation,  which  produces  a  formula  for  the  derivative*  It  is  also  important  to  note 
that  automatic  differentiation  is  not  "numerical  differentiation",  which  yields  only  an 
approximate  value  for  F'(X),  for  example,  by  use  of  the  difference  quotient 

AF  _  F(X  ♦  AX)  -  F(X) 

(3*4j  a* - 5x - ' 

for  a  fixed  value  of  AX. 

As  presented  above,  the  use  of  automatic  differentiation  requires  the  formation  of  the 
additional  code  list  (3.2),  and  evaluation  of  F1 (X)  by  (3.3)  depends  on  values  obtained  in 
the  evaluation  (2.2)  of  F(X).  The  formation  of  an  extra  code  list  is  eliminated,  and  the 
evaluation  of  F(X)  and  F'(X)  are  combined  by  the  use  of  differentiation  arithmetic  applied 
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to  the  original  code  list  (2.1). 

4.  Differentiation  arithmetic.  The  idea  behind  differentiation  arithmetic  it  t  »■ 
familiar  one  of  an  order ed-pair  arithmetic.  Here,  arithmetic  operations  are  defined  in  the 
set  S2  “  S  *  8  of  pairs  of  elements  of  a  set  S  componentwise  in  terms  of  already 
established  operations  in  S.  Some  examples,  with  S  •  K,  the  set  of  integers,  or  S  •  R,  the 
set  of  real  numbers,  are: 

(i)  rational  arithmetic:  X  -  p/q,  p,q  £  K; 

(ii)  complex  arithmetic:  X  «  <x,y)  «  x  +  iy,  x,y  £  R; 

(iii)  interval  arithmetic:  X  »  [a,b]  »  {x  |  a  <  x  <  b,  a,b,x  €  R,  a  <  b). 

The  rules  of  rational  and  complex  arithmetic  are  well-known;  interval  arithmetic  is 
described  in  [1],  [8],  or  [9].  In  any  of  these  arithmetics,  the  function  F(X)  given  by 

(1.1)  can  be  evaluated  by  using  the  code  list  (2.1),  starting  with  a  different  type  of 
element  X,  and  using  different  rules  of  arithmetic  than  in  the  evaluation  of  F(X)  for  a 
real  number  X  «  x  €  R. 

The  basic  elamanta  for  differentiation  arithmetic  are  pairs  of  real  numbers 
F  -  (f,f),  G  -  (g.g’J,  ...  €  D  -  R2. 

For  the  time  being,  the  prism  '  is  to  be  regarded  only  as  a  marker  to  distinguish  the 
second  element  of  a  pair.  The  rules  of  differentiation  arithmetic  are  as  follows: 

1*.  Addition 

(4.1)  F  +  G  -  (f  ,f>  +  (g,g’>  -  (f  ♦  g,  f'  ♦  g'); 

2*.  Subtraction 

(4.2)  F  -  G  -  (f,f)  -  (g,g* )  -  (f  -  g,  f  -  g')j 

3  * .  Multiplication 

(4.3)  F  •  G  -  (f,f')  •  (g,g‘)  “  (f*g,  f*g*  +  g*f')» 

4*.  Division 

(4.4)  F  /  G  -  (f,f)  /  (g,g‘ )  -  (f/g,  g*f'gz  f-*g‘  ),  g  +  0. 


Examination  of  these  rules  shows  their  simple  structure:  The  first  components  of  the 


results  are  formed  by  the  rule  for  evaluation  In  real  arithmetic,  and  the  second  by  the 


corresponding  rule  for  differentiation,  assuming  the  first  components  of  the  operands 
represent  function  values,  and  the  second  components  are  values  of  their  derivatives.  Two 
further  rules  of  differentiation  are  needed  in  order  to  represent  special  quantities  in 
differentiation  arithmetic.  First,  if  x  denotes  a  variable,  the  basic  rule  is 


from  which  it  follows  that  the  independent  variable  X  in  differentiation  arithmetic  is 
represented  by  a  pair 


(4.6) 


X  <*  (x,1). 


where  x  is  the  vales  of  the  independent  variable  X.  Second,  if  c  denotes  a  constant,  then 


(4.7) 


dc 

dx 


0, 


which  means  that  constants  in  differentiation  arithmetic  are  represented  by  pairs 


(4.8) 


C  -  (c,0), 


where  c  is  again  the  value  of  the  constant  C. 

Evaluation  of  a  formula  for  a  rational  function  r  in  differentiation  arithmetic  should 
give  the  correct  values  for  both  F(X)  and  F'(X),  because  this  kind  of  function  can  be 
expressed  as  a  sequence  of  arithmetic  operations  and  the  chain  rule  holds.  To  see  how 
differentiation  arithmetic  works,  try  it  out  on  the  example  (1.1),  following  the  code  list 
(2.1).  Here,  it  is  convenient  to  identify  the  literal  constants  appearing  in  (2.1)  with 
pairs  of  the  form  (4.8).  The  resulting  mixed  arithmetic  (similar  to  multiplication  of 
vectors  by  scalars)  is  less  cumbersome  than  conversion  of  constants  to  the  form  (4.8). 
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With  X  -  (3,1),  (2.1)  becomes 


T1  s-  (3,1)  -  1  «  (2,1)| 

T2  i«  (3,1)  +  3  -  (6,1)1 

(4.9)  T3  :«  (2, 1 ) • (6, 1 )  -  (2*6,  2*1  +  6*1)  «  (12,8)i 
T4  s-  (3,1)  *  2  -  <S,1), 

Ft  (3, 1 ) )  ( 12,8)/(5,  1)  *  (12/5,  (5*8  -  12*1)/52)  =  (12/5,  28/25)i 

so  that 

(4.10)  F(3)  »  12/5  «  2.4,  F'(3)  -  28/25  =  1.12 

for  the  real-valued  function  F(X)  defined  by  (1.1). 

The  rules  of  differentiation  arithmetic  are  actually  simpler  than  the  rules  for 
complex  arithmetic,  and  give  the  correct  values  for  rational  functions  and  their 
derivatives  without  symbolic  manipulations  or  resort  to  numerical  approximations. 
Furthermore,  differentiation  arithmetic  is  easy  to  program  for  hand-held  programmable 
calculators  (see  Appendix  A)  or  computers  (6),  (12),  [13]  (see  Appendix  B) . 

5.  Properties  of  differentiation  arithmetic.  The  reason  why  diffentiation  arithmetic 
works  when  applied  to  rational  functions  is  based  on  the  fact  that  it  has  certain  key 
properties  in  common  with  real  and  complex  arithmetic.  Differentiation  arithmetic  takes 
place  in  the  mathematical  system  D  consisting  of  the  elements  of  R2  with  the  operations 

(4.1) -(4.4).  It  is  easy  to  verify  the  following  properties: 

1*.  Addition  in  D  forms  a  commutative  group.  The  identity  element  for  addition  is 

(5.1)  0  -  (0,0), 

and  each  element  P  €  D  has  the  additive  inverse  (or  negative) 


(5.2) 


-F  -  -tf.f)  -  (-£,  -f ' ) 


2*.  Multiplication  in  D  forms  a  commutative  semigroup  with  an  identity,  that  is, 
multiplication  is  commutative  and  associative,  and 

(5.3)  1  -  (1,0) 

is  the  identity  element  for  multiplication,  F • 1  -  1 *F  -  F  for  all  F  €  O. 

3*.  Multiplication  is  distributive  across  addition, 

(5.4)  F*(G  +  H)  -  F*G  +  F»H, 

for  all  F,G,H  €  D. 

4*.  There  are  no  dlviaors  of  smro,  which  means  that  if 

(5.5)  F*G  »  0, 

and  F  /  0,  then  G  »  0. 

A  mathematical  system  with  these  properties  is  called  an  integral  domain  [4]. 
Furthermore,  every  element  F  -  (f,f)  €  D  such  that  f  0  has  the  Multiplicative  inverse 
(or  reciprocal)  1/F  given  by 

(5.6)  1/F  -  1/(f,f)  -  (1/f,  -f'/f2). 

It  is  ioportant  that  D  is  an  integral  domain  because  this  means  that  the  same  results 
will  be  obtained  independently  of  the  order  in  which  equivalent  sequences  of  arithmetic 
operations  are  performed.  Thus,  for  example,  the  example  function  F(X)  can  be  evaluated  by 
any  code  list  equivalent  to  (2.1)  in  real,  complex,  or  differentiation  arithmetic,  and  the 
results  will  agree.  The  real  number  system  R  and  the  complex  number  system  C  do  share  a 
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property  which  D  does  not  have:  In  R  end  C,  only  the  additive  identity  0  does  not  have  a 
reciprocal;  integral  domains  in  which  this  holds  are  called  fields  [4] . 

6.  functions  in  D.  In  addition  to  rational  functions,  the  student  of  calculus  is 
introduced  to  a  basic  set  of  transcendental  functions.  The  rules  for  differentiation  of 
these  elementary  functions  are  based  on  the  definition  (3.1)  and  its  consequences,  such  as 
the  chain  rule,  i implicit  differentiation,  and  logarithmic  differentiation.  Once  the  rule 
for  differentiation  of  a  function  is  known,  it  can  be  used  to  define  the  corresponding 
function  in  D.  For  exanple, 

(6.1)  eF  »  e^'f'*  ”  (e*,  f*»af), 

(6.2)  In  F  *  In  (f,f)  -  (In  f,  f’/f), 

(6.3)  sin  F  ■  sin  (f ,f • )  -  (sin  f,  f • *coa  f), 

(6.4)  cos  F  -  cos  (f,f)  -  (cos  f,  -  f’*sin  f), 

(6.5)  tan"1  F  «  tan-1  (f,f)  -  (tan_1f,  f '/()  ♦  f2)), 

and  so  on.  The  power  functions  is  usually  handled  in  three  cases: 

(6.6)  Fc  -  <f,f')c  -  (fc,  f'C'f®'1), 

(6.7)  cy  “  c1*'*'1  •  (c*,  f**ln  c  •  cl), 

(6.8)  F®  -  (f,f,)(9,9')  -  (fq,  f * •g»f9"1  +  g*  *ln  g  •  f«), 

inhere  c  €  R  denotes  a  constant. 

By  following  the  pattern  of  the  above  examples,  it  is  easy  to  define  functions  in  D 
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which  give  the  values  of  the  corresponding  function  and  its  derivative  in  R.  In  general, 
by  the  chain  rule, 

(6.9)  u(F)  -  u((f,f))  -  (u(f),  uMfW), 

where  u:  R  ♦  R  is  any  differentiable  function.  The  rules  of  automatic  differentiation  can 
thus  be  extended  to  a  wide  class  of  functions,  and  in  particular  those  differentiable 
functions  encountered  in  elementary  analysis. 

In  Appendix  A,  a  set  of  subroutines  for  the  arithmetic  operations  (4. 1)— (4.4)  and  the 
elementary  functions  (6.1)-(6.5)  are  given  for  a  typical  hand-held  programmable  RPN 
calculator,  the  KP-15C.  A  similar  set  of  operators  and  the  functions  (6. 1)-(6.5)  are  given 
in  Appendix  B  in  the  computer  language  Pascal-SC  [14],  which  can  be  used  on  many 
microcomputers . 

7.  The  geometry  of  differentiation.  It  is  natural  to  identify  the  elements  F  » 

(f,f)  of  D  with  the  points  of  the  real  plane  R2.  In  this  geometry,  constants  lie  on  the 
f-axis  If'  «  0),  and  the  locus  of  the  independent  variable  is  the  line  f'  -  1.  A  function 
P:  D  *  D  can  be  parameterized  by  the  real  variable  x  by  setting  X  “  (x, 1)  to  obtain 

(7.1)  P(X)  »  F(  (X,  1 )  )  -  (f(x),f(x)). 

The  graph  of  F  in  the  f,f -plane  (the  so-called  phase  plane)  reveals  some  interesting  facts 
about  the  underlying  real  function  f:  R  ♦  R.  Points  at  which  the  graph  of  F  crosses  the  f- 
axis  ( f •  «  0)  correspond  to  critical  points  of  the  function  f,  and  crossings  of  the  f’-axis 
(f  ■  0)  to  seros  of  f. 

The  function 

(7.2)  f(x)  -  x2  +  4 

is  represented  by 
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8.  Symbolic  differentiation  revisited.  Symbolic  and  automatic  differentiation  have 
different  goal* •  The  purpose  of  symbolic  differentiation  is  to  produce  a  formula  for  the 
function  f':  R  ♦  K  which  is  the  derivative  of  the  given  function  f:  R  ♦  R,  while  automatic 
differentiation  gives  the  value  f'(x)  of  the  derivative  f'  at  some  given  x  €  R  (as  well  as 
the  value  tlx)).  Both  use  the  same  techniques,  derived  from  the  basic  definition  ( 3. 1 ) • 

It  is  worthwhile  to  note  that  (3,1)  and  the  rules  derived  from  it  actually  define  the 
value  of  the  derivative  F' (X),  and  not  a  formula  for  F ' (X ) .  The  formula  is  obtained  by 
working  with  the  rules  symbolically,  and  using  additional  rules  for  simplification, of  the 
symbolic  results.  With  this  in  mind,  the  rules  of  autosmtic  differentiation  can  also  be 
used  for  symbolic  differentiation.  To  see  this,  evaluate  the  formula  (1,1)  by  using  the 
code  list  (2.1)  with  X  -  (x, 1),  The  results  arei 


(8.1) 


T1  «»  <x,1)  -  1  »  (x  -  1,  1)» 

T2  j«  (x, 1)  +  3  -  (x  +  3,  1)f 

T3  (x-1, 1 ) • (x+3, 1 )  -  <<x-1)*(x+3),  (x-1)*1  +  (x+3)»1) 

-  (x2  +  2x  -  3,  2x  +  2) 

T4  t-  <x,1)  +  2  «  (x  +  2,  1)» 

F(X)  (x2+2x-3,2x+2)/(x+2,1) 

.  (  x2  +  2x  -  3  (x  »  2) »(2x  +  2)  -  (x2  »  2x  -  3)»1  1 

x  +  2  '  (  x  +  2)i  > 

(  x2  *  2x  -  3  x2  +  4x  +  7 
1  x  *  2  '  (x  +  2 ) 2  '  ’ 


This  is  the  right  answer,  obtained  by  applying  the  rules  of  automatic  differentiation 
in  symbolic  form.  This  example  illustrates  the  inherent  difficulty  of  symbolic 
differentiation,  for  example,  in  the  calculation  of  T3,  it  has  to  be  recognized  that  (x- 
1)*1  ♦  (x+3 ) • 1  and  2x  +  2  symbolise  the  same  quantity,  and  the  latter  form  is  preferable. 
In  the  numerical  calculation  of  T3  for  x  ■  3,  the  arithmetic  calculation  2«1  +  6*1  «  8 
required  no  such  simplification.  Computer  programs  which  do  symbolic  differentiation  are 
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large  and  consequently  alow  compared  to  the  small  and  fast  subroutines  required  for 
automatic  differentiation.  This  Is  because  many  extra  rules  for  algebraic  simplification 
must  be  programmed  for  the  symbolic  programs,  and  these  are  not  required  for  simple 
numerical  evaluation  of  functions  and  derivatives. 

9.  Extensions  of  differentiation  arithmetic.  There  are  several  directions  in  which 
the  differentiation  arithmetic  presented  here  for  functions  of  one  real  variable  can  be 
extended.  A  differentiation  arithmetic  can  be  defined  for  higher  derivatives.  In  the  case 
of  second  derivatives,  the  basic  elements  are  triples 

(9.1)  F  -  <f,f',f>, 

the  independent  variable  is  represented  in  the  fora  x  -  (x, 1,0),  and  constants  by  C  * 
(c,0,0).  Leibniz'  rule  is  very  useful  for  the  derivation  of  the  rules  for  multiplication 
and  division  (12).  This  idea  can  be  extended  to  obtain  derivatives  of  arbitrary  order  by 
evaluation  of  the  formula  for  the  original  function,  or  Taylor  coefficients,  since  rules 
are  known  for  formation  of  the  (n+1 ) -tuples 

(9-2)  F  -  (f0,  f, . fn), 

for  arithawtic  operations  and  various  standard  functions,  where 

I*  /.i 

(9.3)  fk  “  kT  f  <x),  k  -  0,1, ...,n, 

with  the  usual  conventions  f^°*(x)  «  f(x)  and  0!  «  1  [12]. 

Another  extension  of  differentiation  arithmetic  is  to  functions  of  several  variables, 
that  is,  to  partial  derivatives.  In  this  case, 

(9.4)  F  -  <f,Vf>, 
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where  Vf  is  the  gradient  vector  of  f,  that  is 


(9.5) 


7f(x) 


f3f(x)  3f(x) 
'  3x^  '  3x 


3f  ( x )  \ 

- 

n 


x  »  (x1 ,x2< • • • .xn) .  The  ith  independent  variable  is  represented  by  “  (x^.e^),  where  e^ 
is  the  ith  unit  vector,  and  constants  by  C  «  (c,0),  where  0  denotes  the  zero  vector  0  * 
(0,0,.. .,0).  The  differentiation  arithmetic  presented  in  this  paper  is  actually  simply 
this  gradient  arithmetic  for  the  special  case  n  »  1. 

Differentiation  arithmetic  can  be  extended  to  higher  partial  derivatives  and  Taylor 
coefficients  in  several  variables  in  the  same  way  as  done  for  functions  of  a  single 
variable.  Second  partial  differentiation  is  based  on  the  representation 


(9.6) 


F  -  <f,Vf,Hf), 


where  Hf  denotes  the  Wsaalan  matrix 

(9.7)  Hf(x,  «  (  |^  ) 

of  second  partial  derivatives  of  f  at  x. 

Increasing  the  order  of  derivatives  needed  or  the  number  of  independent  variables 
quickly  puts  automatic  differentiation  out  of  reach  of  hand  calculation.  However,  the 
required  rules  of  differentiation  arithmetic  are  easy  to  program  for  computers  (5),  [12], 
[13],  and  are  well-suited  to  parallel  computation.  The  results  obtained  have  the  accuracy 
of  symbolic  differention,  but  are  require  much  less  computation.  Furthermore, 
differentiation  arithmetic  can  be  combined  with  other  forms  of  arithmetic,  such  as  interval 
arithmetic,  to  provide  information  about  the  values  of  functions  and  their  derivatives  over 
ranges  of  values  of  the  variables  [8],  [9],  [12],  [13]. 

10.  *  simple  application  of  differentiation  arithmetic.  The  use  of  Newton's  method 

to  solve  an  equation  of  the  form 


(10.1) 


f(x)  -  0 


requires  values  f (xn)  and  f*(xn)  of  the  function  and  its  derivative  for  each  approximation 
xn  of  a  solution  x  •  x*  of  (10.1).  These  values  are  given  directly  by  evaluation  of  the 
function  f  in  differentiation  arithmetic.  For 

(10.2)  X,,  -  (Xn.D, 

the  use  of  differentiation  arithmetic  gives 

(10.3)  Fn  -  F(X„)  -  (f„,  f^>  -  (f  (Xj,),  f(xn)), 

and  thus  for 

(10.4)  x„+1  -  *„  -  fn/fA, 

(10.5)  X„+1  -  (xn+1,1). 

Starting  with  xfl  -  -2.1,  the  results  for  the  function  (1.1)  are: 

xQ  «  2.10000000000 

F(XQ)  -  (27.9000000000,  301.000000000), 
x,  -  -2.19269102990 

F(X,)  -  (13.3762744874,  81.797562426 1), 

X2  -  -2.35622004092 

F(X2)  »  (6.06554082423,  24.6420186899), 
x3  -  -2.60236630436 

F(X3)  -  (2.37799195653,  9.26798946894), 
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2.85894751483 


x4  * 

F(X4)  -  (6.33698948666  x  10”1,  5.06619310633) 
x5  -  -2.98403136682 

F(X5)  -  (6.46519385657  x  10~2,  4.09815663218) 
x6  -  -2.99980722518 

F(X6)  -  (7.71210787890  x  10"4,  4.00115698347) 
x?  -  -2.99999997213 

i 

F(X7)  -  (1.11480002330  x  10“7,  4.00000016722) 

Xg  -  -3.00000000000 

F(Xa)  -  (0.00000000000,  4.00000000000). 

These  results  were  computed  using  the  Paecal-SC  program  given  in  Appendix  C.  This  program 
can  be  modified  easily  to  applied  Newton’s  method  to  functions  other  than  (1.1). 
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APPENDIX  A 


Differentiation  Arithmetic  for  a  Typical  Prograanable  RPN  Calculator 

The  calculator  used  for  these  examples  was  the  HP15-C,  manufactured  by  Hewlett- 
Packard.  The  programing  is  the  same  or  very  similar  for  other  machines  of  this  type. 

The  arithawtic  subroutines  calculate 

(h,h’>  «  (f  ,f '  )*(g,g' ),  where  *  €  {+»-,•./}• 

The  sequence  of  keystrokes  to  enter  F,  G  is  f+f'+g+g',  followed  by  GSB  * ,  where  *  ■  A  for 
addition,  B  for  subtraction,  C  for  multiplication,  and  D  for  division.  The  result  h  will 
be  left  in  X,  and  h'  in  V.  These  routines  use  storage  locations  0,  1  and  2. 


Addition 

Subtraction 

Multiplication 

Division 

f  LBL  A 

f  LBL  B 

f  LBL  C 

f  LBL  D 

X  *  y 

X  ♦  y 

STO  0 

STO  0 

3  TO  0 

STO  0 

R+ 

R+ 

R+ 

R+ 

STO  1 

STO  1 

♦ 

- 

X 

R+ 

x  *  y 

x  *  y 

x  *  y 

x  *  y 

RCL  0 

RCL  0 

STO  x  1 

RCL  1 

+ 

- 

RCL  0 

*■ 

g  RTN 

g  RTN 

X 

♦ 

RCL  1 

g  RTN 

STO  2 

RCL  0 

X 

RCL  1 

* 

RCL  2 

g  RTN 
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The  function  subroutines  calculate 


u<f,f’)  -  (u<f),  u'lfW) 


for  the  functions  u  »  sqrt,  exp.  In,  sin,  cos,  tan-1.  The  sequence  of  keystrokes  is  f  tf ' , 
followed  by  GSB  *,  where  *  -  0  (sqrt),  1  (exp),  .1  (In),  2  (sin),  3  (cos),  4  (tan-1).  The 
result  u(f)  will  be  left  in  X,  and  the  derivative  u* (f)*f*  will  be  left  in  Y.  These 
routines  use  storage  location  0. 


aart(F)  -  /F 

sxd(F)  -  eF 

ln(F) 

f  LBL  0 

f  LBL  1 

f  LBL 

X  ♦  Y 

X  *  Y 

X  *  Y 

/x 

,  •* 

STO  0 

STO  0 

STO  0 

♦ 

♦ 

X 

RCL  0 

2 

RCL  0 

g  LH 

♦ 

g  RTN 

g  RTH 

RCL  0 

g  BTH 

»L"<F2 

£°ELEL 

tan”,(F) 

f  LBL  2 

f  LBL  3 

f  LBL  4 

g  RAO 

g  RAO 

g  RAO 

X  *  Y 

X  *  Y 

X  *  Y 

STO  0 

STO  0 

STO  0 

COS 

SIN 

g  x2 

X 

X 

1 

RCL  0 

CHS 

♦ 

SIN 

RCL  0 

* 

g  RTH 

COS 

RCL  0 

g  RTH 

g  tan"1 
g  RTH 

-21- 


APPENDIX  B 


Pascal-SC  Operators  and  Functions  for  Differentiation  Arithmetic 


Differentiation  arithmetic  is  performed  in  Pascal-SC  by  use  of  the  data  type 
DERIV,  which  is  declared  as  follows: 

TYPE  DERIV  -  RECORD  X, PRIME:  REAL  END; 

This  means  that  if  F  is  a  variable  of  type  DERIV,  then  F.X  =  F(X)  represents  its 
function  value,  and  F. PRIME  -  F‘(X)  the  value  of  its  derivative.  The  independent 
variable  X  has  X. PRIME  »  1,  while  C. PRIME  =  0  for  constants  C. 

In  order  to  be  able  to  mix  REAL,  INTEGER,  and  DERIV  types  in  expressions,  where 
REAL  and  INTEGER  are  considered  to  be  constants,  22  arithmetic  operators  are  r 
required.  If  K,  R,  and  D  denote  respectively  generic  variables  of  types  INTEGER, 
REAL,  and  DERIV,  then  the  following  operators  are  needed: 

Addition:  +D,  K  +  D,  D  +  K,  R  +  D,  D  +  R,  D  +  D; 

Subtraction:  -D,  K  -  D,  D  -  K,  R  -  D,  D  -  R,  D  -  D; 

Multiplication:  K*D,  D*K,  R*D,  D*R,  D*D; 

Division:  K/D,  D/K,  R/D,  D/R,  D/D. 

The  following  files  contain  the  source  code  for  the  above  operators: 

Addition  (DADD. OPR),  subtraction  (DSUB.OPR),  multiplication  (DMUL.OPR),  d 
division  (DDIV.OPR) . 


1.  ADDITION  OPERATORS 

OPERATOR  +  (A:  DERIV)  RES:  DERIV; 

BEGIN 
RES  :-  A 
END; 

OPERATOR  +  (K:  INTEGER; A:  DERIV)  RES:  DERIV; 
BEGIN 

A.X  :-  K  +  A.X; 

RES  :-  A 
END; 

OPERATOR  +  (A:  DERIV ;K:  INTEGER)  RES:  DERIV; 
BEGIN 

A.X  :-  A.X  +  K; 

RES  :-  A 
END; 
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OPERATOR  +  (R:  REAL ; A :  DERIV)  RES:  DERIVj 
BEGIN 

A.X  :-  R  +  A.X; 

RES  :=  A 
END; 

OPERATOR  +  (A:  DERIV ;R: REAL)  RES:  DERIV; 
BEGIN 

A.X  :=  A.X  +  R; 

RES  :  =  A 
END; 

OPERATOR  +  (A,B;  DERIV)  RES:  DERIV; 

BEGIN 

A.X  :«  A.X  +  B.X; 

A. PRIME  :«  A. PRIME  +  B. PRIME; 

RES  :-  A; 

END; 


2 .  SUBTRACTION  OPERATORS 

OPERATOR  -  (A:  DERIV)  RES:  DERIV; 

BEGIN 

A.X  :-  -A.X; 

A. PRIME  :-  -A. PRIME; 

RES  :-  A 
END; 

OPERATOR  -  (K:  INTEGER; A:  DERIV)  RES:  DERIV; 
BEGIN 

A.X  :«  K  -  A.X) 

A. PRIME  :-  -A. PRIME; 

RES  :-  A 
END; 

OPERATOR  -  (A:  DERIV ;K:  INTEGER)  RES:  DERIV; 
BEGIN 

A.X  :-  A.X  -  K; 

RES  :-  A 
END; 

OPERATOR  -  (R:  REAL; A:  DERIV)  RES:  DERIV; 
BEGIN 

A.X  :-  R  -  A.X; 

A. PRIME  :-  -A. PRIME; 

RES  :»  A 
END; 

OPERATOR  -  (A:  DERIV ;R:  REAL)  RES:  DERIV; 
BEGIN 

A.X  :-  A.X  -  R; 

RES  :-  A 
END; 
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OPERATOR  -  (A,B:  DERIV)  RES:  DERIVj 
BEGIN 

A.X  :«  A.X  -  B.X; 

A. PRIME  A. PRIME  -  B. PRIME; 
RES  :•  A 
END; 


3 .  MULTIPLICATION  OPERATORS 

OPERATOR  *  <K:  INTEGER; A:  DERIV)  RES:  DERIV; 
BEGIN 

A.X  :-  K*A.X; 

A. PRIME  K*A. PRIME; 

RES  :  =  A; 

END; 

OPERATOR  *  (A:  DERIV ;K:  INTEGER)  RES:  DERIV; 
BEGIN 

A.X  :«  A . X*K ; 

A. PRIME  :»  A.PRIME*K; 

RES  :>  A; 

END; 

OPERATOR  *  (R:  REAL; A:  DERIV)  RES:  DERIV; 
BEGIN 

A.X  :»  R*A.X; 

A. PRIME  .•=  R*A. PRIME; 

RES  :•  A; 

END; 

OPERATOR  *  (A:  DERIV; R:  REAL)  RES:  DERIV; 
BEGIN 

A.X  :-  A . X*R ; 

A. PRIME  A . PRIME*R ; 

RES  A 
END; 

OPERATOR  *  (A,B:  DERIV)  RES:  DERIV; 

VAR  U:  DERIV; 

BEGIN 

U.X  A . X*B .X; 

U. PRIME  :-  A.X*B.PRIME  +  A.PRIME*B.X; 
RES  U 
END; 


4.  DIVISION  OPERATORS 

OPERATOR  /  (K:  INTEGER; A:  DERIV)  RES:  DERIV; 
VAR  U:  DERIV; 

BEGIN 

U.X  K/A.X; 

U. PRIME  - A . PRIME*U . X/ A . X ; 

RES  :-  U 
END; 
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OPERATOR  /  (A:  OERIVjK:  INTECER)  RE S:  DERIV; 
BEGIN 

A.X  A.X/K; 

A. PRIME  A.PRIME/K; 

RES  A 
END; 

OPERATOR  /  (R:  REAL; A:  DERIV)  RES:  DERIV; 

VAR  U:  DERIV; 

BEGIN 

U.X  R/A.X; 

U. PRIME  -  A .  PRIME*!! .  X/ A .  X ; 

RES  :•  U; 

END; 

OPERATOR  /  (A:  DERIV; R:  REAL)  RES:  DERIV; 
BEGIN 

A.X  :-  A.X/R; 

A. PRIME  :-  A.PRIME/R; 

RES  A 
END; 

OPERATOR  /  (A,B:  DERIV)  RES:  DERIV; 

VAR  V:  DERIV; 

BEGIN 

U.X  A.X/B.X; 

U. PRIME  (A. PRIME  -  U . X+B .PRIME) /B .X; 
RES  :-  U 
END; 


Functions  for  type  DERIV  have  names  beginning  with  D  and  source  code  for  them  is 
given  in  the  files  *.FUN,  where  *  »  DSQRT,  DEXP,  DLN,  DSIN,  DCOS,  DARCTAN. 

1.  SQUARE  ROOT 

FUNCTION  DSQRT(A:  DERIV):  DERIV; 

BEGIN 

A.X  SQRT(A.X); 

A. PRIME  :«  0.5*A. PRIME/A. X; 

DSQRT  :«  A 
END; 


2.  EXPONENTIAL  FUNCTION  (base  e) 

FUNCTION  DEXP(A:  DERIV):  DERIV; 
BEGIN 

A.X  :-  EXP(A.X) ; 

A. PRIME  A.X*A.PRIME; 

DEXP  A 
END; 
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3 .  NATURAL  LOGARITHM 

FUNCTION  DLN(A:  DERIV) :  DERIV; 
VAR  U:  DERIV; 

BEGIN 

U.X  i-  LN(A.X); 

U. PRIME  !-  A.PRIME/A.Xj 
DLN  :-  U 
END; 


A.  SINE 

FUNCTION  DSIN(A:  DERIV):  DERIV; 
VAR  U:  DERIV; 

BEGIN 

U.X  SIN(A.X); 

U. PRIME  A . PRIME*COS ( A . X ) ; 

DSIN  :-  U 
END; 


S.  COSINE 

FUNCTION  DCOS(A:  DERIV):  DERIV; 
VAR  U:  DERIV; 

BEGIN 

U. X  :-  COS(A.X); 

V.  PRIME  .—  -A.PRINE*SIN(A.X>; 
DCOS  :■  u 

END; 


6 .  ARCTANGENT 


FUNCTION  DARCTAN(A:  DERIV):  DERIV; 

VAR  U:  DERIV; 

BEGIN 

U.X  :-  ARCTAN(A.X) ; 

U. PRIME  :-  A. PRIME/ (1.0  +  A.X*A.X) 
DARCTAN  :-  U 
END; 


APPENDIX  C 


A  Pascal 'SC  Program  for  Newton's  Method 


PROGRAM  LITTLEN(INPUT, OUTPUT); 

(*  This  program  applies  Newton's  method  to  find  a  zero  of  a  function 
f(x)  of  one  variable  using  differentiation  arithmetic.  *) 

TYPE  DERIV  -  RECORD  X, PRIME;  REAL  END; 

(*  If  F  Is  of  type  DERIV,  then  F.X  -  F(X),  F. PRIME  -  F'(X).  *) 

VAR  X:  DERIV;  (*  THE  INDEPENDENT  VARIABLE  *) 

F:  DERIV;  (*  THE  DEPENDENT  VARIABLE  *) 

K:  INTEGER;  (*  ITERATION  COUNTER  *) 

C:  CHAR;  (*  CONTROL  CHARACTER  *) 

(*  INSERT  THE  OPERATORS  AND  FUNCTIONS  FOR  TYPE  DERIV  YOU  NEED  TO 
EVALUATE  YOUR  FUNCTION  HERE,  FOR  EXAMPLE:  *) 

$ INCLUDE  DADD. OPR; 

| INCLUDE  DSUB.OPR; 

$ INCLUDE  DMUL.OPR; 

♦INCLUDE  DDIV.OPR; 

(*  END  OF  DERIV  OPERATORS  AND  FUNCTIONS;  THE  ONES  GIVEN  ABOVE  ARE 
ADEQUATE  FOR  ANY  RATIONAL  FUNCTION.  *) 

BEGIN  (*  PROGRAM  LITTLEN  *) 

X. PRIME  :«  1.0;C  'R';WHILE  C  -  'R'  DO 

BEGIN  (*  SET  INITIAL  VALUES  *) 

K  0;  (*  SET  ITERATION  COUNTER  TO  0  *) 

WRITELN ; 

WRITELN( 'Enter  Initial  value  of  X; ’ ) ;READ(X.X) ; 

C  ' I ' ; WHILE  C  -  'I*  DO 

BEGIN  (*  NEWTON  ITERATION  *} 

(*  INSERT  A  STATEMENT  OR  STATEMENTS  TO  EVALUATE  F(X)  HERE, 

FOR  EXAMPLE;  A) 

F  (X  -  1)A(X  +  3)/(X  +  2); 


(*  END  OF  DEFINITION  OF  F(X)  *) 


SECURITY  CLASSIFICATION  OF  THIS  RAOK  fWhwi  Pi*  inland) 


REPORT  DOCUMENTATION  PAGE 

READ  INSTRUCTIONS 

BEFORE  COMPLETING  FORM 

1.  REFORT  NUMBER  2.  OOVT  ACCESSION  NO. 

2688 

S.  RECIPIENT'S  CATALOG  NUMBER 

4.  TITLE  (and  Subtltt*) 

THE  ARITHMETIC  OP  DIFFERENTIATION 

S.  TYPE  OF  REPORT  •  PERlOO  COVERED 

Summary  Report  -  no  specific 
reporting  period 

S.  PERFORMING  ORG.  REPORT  NUMBER 

7-  author!*; 

L.  B.  Rail 

4.  CONTRACT  OR  GRANT  NUMBER!*! 

DAAG29-80-C-00  41 

S.  FERFORMINO  ORGANIZATION  NAME  AND  ADDRESS 

Mathematics  Research  Center,  University  of 

610  Walnut  Street  Wisconsin 

Madison.  Wisconsin  53706 

l6.  PROGRAM  ELEMENT.  PROJECT,  TASK 

AREA  k  WORK  UNIT  NUMBERS 

Work  Unit  Number  3  - 
Numerical  Analysis  and 
Scientific  Computing 

II.  CONTROLLING  OFFICE  NAME  AND  ADDRESS 

U.  S.  Army  Research  Office 

P.O.  Box  12211 

Research  Triangle  Park.  North  Carolina  27709 

12.  REPORT  DATE 

May  1984 

IS.  NUMBER  OF  PAGES 

28 

14.  MONITORING  AGENCY  name  A  ADDRESS!!!  dlMannl  boat  Control  1  Inf  Ottleo) 

IS.  SECURITY  CLASS,  (al  till*  man) 

UNCLASSIFIED 

IS-  DISTRIBUTION  STATEMENT  (al  ml*  R*patl) 

Approved  for  public  release;  distribution  unlimited. 

•7.  DISTRIBUTION  STATEMENT  (at  d>*  abotracl  inland  hi  Block  20.  II  efferent  1 root  Rapart) 

IS.  SUPPLEMENTARY  NOTES 

IS.  KEY  WORDS  (Canllnua  an  *******  aid*  It  nacaaaar*  mid  Idtntlly  black  numb**) 

Automatic  differentiation 

20.  ABSTRACT  fC enttnue  on  *******  old*  II  neceeeeqr  and  Identify  br  block  numb**) 

This  report  describes  automatic  differentiation,  which  is  neither 

symbolic  nor  approximate,  for  single  functions  of  one  real  variable.  The 

rules  of  evaluation  and  differentiation  are  combined  into  an  ordered-pair 

arithmetic  similar  to  complex  arithmetic,  but  slightly  simpler.  Evaluation 

of  the  formula  for  a  function  in  this  arithmetic  yields  both  the  values  of 

(cont. ) 

DO  ,  jSTn  1473  «oition  of  I  nov  m  is  oosolctk  UNCLASSIFIED 


SECURITY  CLASSIFICATION  OF  THIS  FAOE  fWMn  Dm*  inland) 

I 


ABSTRACT  (cont.) 


the  function  and  its  derivative,  without  a  formula  for  the  derivative  of  the 
function,  and  without  numerical  approximation,  since  this  arithmetic  is  based 
on  the  well-known  rules  for  differentiation.  The  properties  of  this 
arithmetic  are  examined,  and  illustrated  by  simple  examples.  Subroutines 
are  given  for  differentiation  arithmetic  both  on  a  hand-held  programmable 
calculator,  and  in  the  microcomputer  language  Pascal-SC.  An  application  of 
this  arithmetic  to  the  solution  of  equations  by  Newton's  method  is  given, 
using  a  Pascal-SC  program. 


